# Title     : TODO
# Objective : TODO
# Created by: Administrator
# Created on: 2019/8/12
library(optparse)
library(magrittr)
library(impute)
library(tidyverse)
library(NormalizeMets)
library(FitAR)

option_list <- list(
  make_option("--cc", default = "calculate_config.txt", type = "character", help = "config file"),
  make_option("--si", default = "1;0", type = "character", help = "step and config arg index"),
  make_option("--g", default = "group.txt", type = "character", help = "sample group file")
)
opt <- parse_args(OptionParser(option_list = option_list))

library(limma)
library(dplyr)

# Create the raw data matrix.
set.seed(16)
mat <- matrix(rnorm(900), 100, 9)
colnames(mat) <- str_c('Sample', 1:9)
rownames(mat) <- paste('Var', 1:nrow(mat))
matOut<-mat %>%
  as.data.frame() %>%
  rownames_to_column("Metabolite")
write_csv(matOut, "mat.csv")

## Reorder the Sample inform data.
age <- c(45, 55, 47, 56, 43, 56, 66, 68, 69)
gender <- factor(c("M", "F", "F", "F", "M", "F", "F", "M", "F"))
condition <- factor(c(rep("healthy", 3), rep("pre_cancer", 3), rep("cancer", 3)),
                    levels = c("healthy", "pre_cancer", 'cancer'))
sample_info_df <- data.frame(Sample_Info = str_c('Sample', 9:1), age, gender, condition)
sample_info_df <- match(colnames(mat), sample_info_df$Sample_Info) %>% sample_info_df[.,]
sample_info_df
write_csv(sample_info_df, "sample_info.csv")

contrast_note1 <- "conditionpre_cancer-conditionhealthy"

design_D <- model.matrix(data = sample_info_df, ~0 + condition + age + gender)
print("=in=")
design_D
print("=out=")

cont.matrix_D <- makeContrasts(contrast1 =
                                 # "conditionpre_cancer-conditionhealthy",
                                 conditionpre_cancer - conditionhealthy,
                               # "conditionpre_cancer-conditionhealthy",
                               contrast2 =
                                 'conditioncancer-conditionhealthy',
                               levels = design_D)
cont.matrix_D

cont.matrix_D %>% class

fit_D <- lmFit(mat, design_D)
fit_D <- contrasts.fit(fit_D, cont.matrix_D)
fit_D <- eBayes(fit_D)
fit_D

Treat1.adj <- topTable(fit_D, coef = "contrast1", number = nrow(mat), adjust.method = 'BH')
Treat2.adj <- topTable(fit_D, coef = "contrast2", number = nrow(mat), adjust.method = 'BH')
# 去除log2FC列。
Treat2.adj <- Treat2.adj[, colnames(Treat2.adj) != 'logFC']

tt_D <- topTable(fit_D, coef = 0)
print("=in=")

Treat1.adj <- Treat1.adj[, colnames(Treat1.adj) != 'logFC']
Treat1.adj %>% head
Treat1.adj <- Treat1.adj[, c('t', 'P.Value', 'adj.P.Val')]
colnames(Treat1.adj) <- c('t', "P.value", 'FDR')
Treat1.adj %>% write.csv(., 'demo res.csv')
mean(mat[1,])
